############################################
# Replication Script for Religious Roots of Belief in Misinformation
# Authors: Sumitra Badrinathan and Simon Chauchard
# Date: May 21, 2025
############################################

########### SECTION 1: SET UP ###########
# install packages 

packages <- c("dotwhisker", "effects", "ggeffects", "stargazer", 
              "margins", "ggplot2", "tidyverse", "here")

install_if_missing <- function(pkg) {
  if (!requireNamespace(pkg, quietly = TRUE)) install.packages(pkg, dependencies = TRUE)
  library(pkg, character.only = TRUE)
}
lapply(packages, install_if_missing)

# load data 
covid <- read.csv("covid-misinformation.csv")


########## SECTION 2: CREATE VARIABLES  #######
##### conspiracy misinfo DV variable #####

# recoding such that 1 = correct discernment
table(covid$B6_Q8)
covid$jamaat <- 0
covid$jamaat[covid$B6_Q8==3 | covid$B6_Q8==4] <- 1
table(covid$jamaat)

table(covid$B6_Q9)
covid$sneezing <- 0
covid$sneezing[covid$B6_Q9==3 | covid$B6_Q9==4] <- 1
table(covid$sneezing)

table(covid$B6_Q10)
covid$bio <- 0
covid$bio[covid$B6_Q10==3 | covid$B6_Q10==4] <- 1
table(covid$bio)

table(covid$B6_Q11)
covid$foreign <- 0
covid$foreign[covid$B6_Q11==3 | covid$B6_Q11==4] <- 1
table(covid$foreign)

table(covid$B6_Q12)
covid$lab <- 0
covid$lab[covid$B6_Q12==1 | covid$B6_Q12==2] <- 1
table(covid$lab)

table(covid$B6_Q13)
covid$fiveg <- 0
covid$fiveg[covid$B6_Q13==1 | covid$B6_Q13==2] <- 1
table(covid$fiveg)

# adding them up
covid$ConspiracyMisinfo <- covid$jamaat + covid$sneezing + covid$bio + covid$foreign +
  covid$lab + covid$fiveg

##### medical misinfo DV variable #####

table(covid$B5_Q8)
covid$homeo <- 0
covid$homeo[covid$B5_Q8==3 | covid$B5_Q8==4] <- 1
table(covid$homeo)

table(covid$B5_Q9)
covid$kalonji <- 0
covid$kalonji[covid$B5_Q9==3 | covid$B5_Q9==4] <- 1
table(covid$kalonji)

table(covid$B5_Q10)
covid$genetic <- 0
covid$genetic[covid$B5_Q11==3 | covid$B5_Q11==4] <- 1
table(covid$genetic)

table(covid$B5_Q11)
covid$breath <- 0
covid$breath[covid$B5_Q11==3 | covid$B5_Q11==4] <- 1
table(covid$breath)

table(covid$B5_Q12)
covid$antib <- 0
covid$antib[covid$B5_Q12==1 | covid$B5_Q12==2] <- 1
table(covid$antib)

table(covid$B5_Q13)
covid$bleach <- 0
covid$bleach[covid$B5_Q13==1 | covid$B5_Q13==2] <- 1
table(covid$bleach)

# adding them up
covid$CuresMisinfo <- covid$homeo + covid$kalonji + covid$genetic + covid$breath +
  covid$antib + covid$bleach

##### combined DV for all stories  #####

covid$AllMisinfo <- covid$CuresMisinfo + covid$ConspiracyMisinfo
summary(covid$AllMisinfo)

##### religiosity scale #####
# recoding such that higher = more religious

# would marry someone who is not hindu, 1=strongly agree
table(covid$B3_Q1) # already coded correctly

# an athiest can be moral, 1=strongly agree
table(covid$B3_Q8) # already coded correctly

# religion can help me cope better
table(covid$B3_Q2)
covid$relig1 <- NA
covid$relig1[covid$B3_Q2==4] <- 1
covid$relig1[covid$B3_Q2==3] <- 2
covid$relig1[covid$B3_Q2==2] <- 3
covid$relig1[covid$B3_Q2==1] <- 4
table(covid$relig1)

# important to teach children
table(covid$B3_Q3)
covid$relig2 <- NA
covid$relig2[covid$B3_Q3==4] <- 1
covid$relig2[covid$B3_Q3==3] <- 2
covid$relig2[covid$B3_Q3==2] <- 3
covid$relig2[covid$B3_Q3==1] <- 4
table(covid$relig2)

# fasting
table(covid$B3_Q4)
covid$relig3 <- NA
covid$relig3[covid$B3_Q4==4] <- 1
covid$relig3[covid$B3_Q4==3] <- 2
covid$relig3[covid$B3_Q4==2] <- 3
covid$relig3[covid$B3_Q4==1] <- 4
table(covid$relig3)

# god blesses me puja
table(covid$B3_Q5)
covid$relig4 <- NA
covid$relig4[covid$B3_Q5==4] <- 1
covid$relig4[covid$B3_Q5==3] <- 2
covid$relig4[covid$B3_Q5==2] <- 3
covid$relig4[covid$B3_Q5==1] <- 4
table(covid$relig4)

# consult w astrologer
table(covid$B3_Q6)
covid$relig5 <- NA
covid$relig5[covid$B3_Q6==4] <- 1
covid$relig5[covid$B3_Q6==3] <- 2
covid$relig5[covid$B3_Q6==2] <- 3
covid$relig5[covid$B3_Q6==1] <- 4
table(covid$relig5)

# vegetarian food
table(covid$B3_Q7)
covid$relig6 <- NA
covid$relig6[covid$B3_Q7==4] <- 1
covid$relig6[covid$B3_Q7==3] <- 2
covid$relig6[covid$B3_Q7==2] <- 3
covid$relig6[covid$B3_Q7==1] <- 4
table(covid$relig6)

# adding
covid$religion_added <- covid$B3_Q1 + covid$B3_Q8 + covid$relig1 + covid$relig2 +
  covid$relig3 + covid$relig4 + covid$relig5 + covid$relig6 
summary(covid$religion_added)
covid <- covid[!is.na(covid$religion_added), ] 

# standardized scale
covid$religiosity <- (covid$religion_added - mean(covid$religion_added))/ sd(covid$religion_added)
summary(covid$religiosity)


##### covariates #####

# gender
covid$gender[covid$B1_Q3==1] <- "Male"
covid$gender[covid$B1_Q3==2] <- "Female"
table(covid$gender)

# indicator for upper caste
covid$uppercaste <- 0
covid$uppercaste[covid$B1_Q7==1] <- 1
table(covid$uppercaste)

# education
covid$education <- NA
covid$education[covid$B1_Q4<5 | covid$B1_Q4==9] <- 1
covid$education[covid$B1_Q4==5 | covid$B1_Q4==6] <- 2
covid$education[covid$B1_Q4==7 | covid$B1_Q4==8] <- 3
table(covid$education)

# bjp support binarized
covid$BJP <- 0
covid$BJP[covid$B3_Q14==1 | covid$B3_Q14==2] <- 1
table(covid$BJP)

# science knowledge scale
covid$science1 <- 0 
covid$science1[covid$B4_Q1==2] <- 1

covid$science2 <- 0 
covid$science2[covid$B4_Q2==1] <- 1

covid$science3 <- 0 
covid$science3[covid$B4_Q3==1] <- 1

covid$science4 <- 0 
covid$science4[covid$B4_Q4==1] <- 1

covid$science5 <- 0 
covid$science5[covid$B4_Q5==1] <- 1

covid$science6 <- 0 
covid$science6[covid$B4_Q6==1] <- 1

covid$science7 <- 0 
covid$science7[covid$B4_Q7==1] <- 1

covid$science8 <- 0 
covid$science8[covid$B4_Q8==2] <- 1

covid$science_knowledge <- covid$science1 + covid$science2 + covid$science3 + covid$science4 +
  covid$science5 + covid$science6 + covid$science7 + covid$science8


############## SECTION 3: MAIN TEXT RESULTS ###########

##### figure 2 #####
conspiracies <- rbind(data.frame(prop.table(table(covid$jamaat))), 
                      data.frame(prop.table(table(covid$sneezing))),
                      data.frame(prop.table(table(covid$bio))),
                      data.frame(prop.table(table(covid$foreign))),
                      data.frame(prop.table(table(covid$lab))),
                      data.frame(prop.table(table(covid$fiveg))))
#subsetting to those who did not correctly discern
conspiracies <- subset(conspiracies, conspiracies$Var1=="0")
conspiracies$story <- c("Jamaat a conspiracy\n to spread covid", "Muslim devotees \n sneezing to spread covid", "Covid a Chinese\nbiowarfare weapon", 
                        "Foreign powers\nspread covid", 
                        "Covid not\n created in a lab", "Covid not \n linked to 5G")
conspiracies$Freq <- conspiracies$Freq*100
conspiracies$Freq <- round(conspiracies$Freq, digits=2)
conspiracies$type <- rep("Conspiracy Theory", nrow(conspiracies))
conspiracies$veracity <- c("False Headlines", "False Headlines", "False Headlines", "False Headlines", "True Headlines", "True Headlines")

#  cures
cures <- rbind(data.frame(prop.table(table(covid$homeo))), 
               data.frame(prop.table(table(covid$kalonji))),
               data.frame(prop.table(table(covid$genetic))),
               data.frame(prop.table(table(covid$breath))),
               data.frame(prop.table(table(covid$antib))),
               data.frame(prop.table(table(covid$bleach))))
cures <- subset(cures, cures$Var1=="0")
cures$story <- c("Homeopathy\nprevents covid", "Kalonji seeds\ncure covid", "Indians genetically\nprotected from covid", 
                 "Test for covid\nby holding breath", "Antibiotics do not\ncure covid", 
                 "Bleach does not\nprotect from covid")
cures$Freq <- cures$Freq*100
cures$Freq <- round(cures$Freq, digits=2)
cures$type <- rep("Miracle Cure", nrow(cures))
cures$veracity <- c("False Headlines", "False Headlines", "False Headlines", "False Headlines", "True Headlines", "True Headlines")

allbars <- rbind(cures, conspiracies)
allbars$Freq <- round(allbars$Freq)

q <- ggplot(data=allbars, aes(x=reorder(story, Freq), y=Freq, fill=veracity)) +
  geom_bar(stat = "identity", width=0.8) +
  geom_text(aes(label = Freq), hjust=-0.2, vjust=0.5, size=3.4) 
q  + coord_flip() + theme_bw() +
  xlab("") + ylab("Percent Incorrectly Classifying Headline") + 
  ggtitle("") +
  scale_fill_grey(start=0.4, end=0.7) +
  theme(legend.position="bottom", legend.title = element_blank())

##### figure 3 #####
hyp1b <- lm(AllMisinfo ~ religiosity, data=covid)
summary(hyp1b)

cplot(hyp1b, "religiosity", what = "prediction", 
      main = "", 
      xlab="Religiosity (Standardized Scale)", 
      ylab="Predicted # of stories correctly identified")


##### table 1 #####
covid$B5_Q1 <- factor(covid$B5_Q1, levels = c("5", "1", "2", "3", "4"))
levels(covid$B5_Q1)

main_cures <- lm(CuresMisinfo ~ factor(B5_Q1), data=covid)

main_cons <- lm(ConspiracyMisinfo ~ factor(B5_Q1), data=covid)

stargazer(main_cons, main_cures, star.cutoffs =  c(0.05, .001, 0.001), type="text")

##### table 2 #####
# CONSPIRACY

# conspiracy true (recode as higher = more accurate)
table(covid$B6_Q12) # lab
table(covid$B6_Q13) # 5g
covid$labnew <- 5-covid$B6_Q12 # higher (4) = more accurate
covid$fivegnew <- 5-covid$B6_Q13

# mean of true stories 
covid$cons_true_mean <- (covid$labnew + covid$fivegnew) / 2 

# z scores true stories
covid$cons_true_z <- (covid$cons_true_mean - mean(covid$cons_true_mean)) / sd(covid$cons_true_mean)

# conspiracy false stories mean
covid$B6_Q8new <- 5-covid$B6_Q8
covid$B6_Q9new <- 5-covid$B6_Q9
covid$B6_Q10new <- 5-covid$B6_Q10
covid$B6_Q11new <- 5-covid$B6_Q11

covid$cons_false_mean <- (covid$B6_Q8new + covid$B6_Q9new + covid$B6_Q10new + covid$B6_Q11new) / 4

# z scores of mean of false stories
avg <- mean(covid$cons_false_mean, na.rm=T)
stdev <-  sd(covid$cons_false_mean, na.rm=T)
covid$cons_false_z <- (covid$cons_false_mean - avg) / stdev

# discernment conspiracies
covid$discern_cons_z <- covid$cons_true_z - covid$cons_false_z

# CURES

table(covid$B5_Q12) # antibiotics
table(covid$B5_Q13) # bleach
covid$antibnew <- 5-covid$B5_Q12 # higher (4) = more accurate
covid$bleachnew <- 5-covid$B5_Q13

# mean of true stories 
covid$cures_true_mean <- (covid$antibnew + covid$bleachnew) / 2 

# z scores true stories
covid$cures_true_z <- (covid$cures_true_mean - mean(covid$cures_true_mean)) / sd(covid$cures_true_mean)

# cures false stories mean
covid$B5_Q8new <- 5- covid$B5_Q8
covid$B5_Q9new <- 5- covid$B5_Q9
covid$B5_Q10new <- 5- covid$B5_Q10
covid$B5_Q11new <- 5- covid$B5_Q11

covid$cures_false_mean <- (covid$B5_Q8new + covid$B5_Q9new + covid$B5_Q10new + covid$B5_Q11new) / 4

# z scores false stories
avg2 <- mean(covid$cures_false_mean, na.rm=T)
stdev2 <-  sd(covid$cures_false_mean, na.rm=T)
covid$cures_false_z <- (covid$cures_false_mean - avg2) / stdev2

# discernment cures
covid$discern_cures_z <- covid$cures_true_z - covid$cures_false_z

#models
covid$B5_Q1 <- factor(covid$B5_Q1, levels = c("5", "1", "2", "3", "4"))

regdiscerncons <- lm(discern_cons_z ~ factor(B5_Q1), data=covid)
regdiscerncure <- lm(discern_cures_z ~ factor(B5_Q1), data=covid)

stargazer(regdiscerncons, regdiscerncure, star.cutoffs = c(0.05, 0.01, 0.001), type="text")

##### table 3 #####
covid$B5_Q1 <- factor(covid$B5_Q1, levels = c("4", "1", "2", "3", "5"))
levels(covid$B5_Q1)

main_cures2 <- lm(CuresMisinfo ~ factor(B5_Q1), data=covid)

main_cons2 <- lm(ConspiracyMisinfo ~ factor(B5_Q1), data=covid)

stargazer(main_cons2, main_cures2, star.cutoffs =  c(0.05, .001, 0.001), type="text")


# APPENDIX




############## SECTION 4: APPENDIX RESULTS ############

##### table D1 #####
stats <- data.frame(covid$B5_Q1, covid$religiosity, covid$BJP, covid$B1_Q3,
                    covid$B1_Q2, covid$gender, covid$B1_Q8, covid$education, covid$uppercaste, 
                    covid$B2_Q7, covid$B2_Q8, covid$science_knowledge)

stargazer(stats, summary.stat = c("n", "mean", "sd", 
                                  "min", "median", "max"), type="text")


##### table E1 #####
covid$B5_Q1 <- factor(covid$B5_Q1, levels = c("5", "1", "2", "3", "4"))
levels(covid$B5_Q1)
cures_controls <- lm(CuresMisinfo ~ factor(B5_Q1) + religiosity + BJP + science_knowledge +
                       B1_Q2 + gender + B1_Q8 + education + uppercaste + B2_Q7 + 
                       B2_Q8, data=covid)
summary(cures_controls)

consp_controls <- lm(ConspiracyMisinfo ~ factor(B5_Q1) + religiosity + BJP +
                       B1_Q2 + gender + B1_Q8 + education + uppercaste +
                       science_knowledge + B2_Q7 + B2_Q8, data=covid)
summary(consp_controls)

stargazer(consp_controls, cures_controls, star.cutoffs =  c(0.05, .001, 0.001), type="text")


##### table F1 #####

covid$check1 <- 0
covid$check1[covid$fltrr2==1 & covid$fltrr4==1] <- 1
table(covid$check1)

covid$check2 <- 0
covid$check2[covid$B2_Q11==3] <- 1
table(covid$check2)

# both checks combined
covid$both_checks <- covid$check1 + covid$check2

covid$B5_Q1 <- as.factor(covid$B5_Q1)
covid$B5_Q1 <- relevel(covid$B5_Q1, "5")

# reproducing main effect w checks variable as continuous
main_cures_attn <- lm(CuresMisinfo ~ factor(B5_Q1) + both_checks, data=covid)

main_cons_attn <- lm(ConspiracyMisinfo ~ factor(B5_Q1) +  both_checks, data=covid)

stargazer(main_cons_attn, main_cures_attn, star.cutoffs =  c(0.05, .001, 0.001), type="text")


##### table G1 #####
hyp2 <- subset(covid, covid$B5_Q1==1 | covid$B5_Q1==5)
hyp2$treatment <- 0
hyp2$treatment[hyp2$B5_Q1==1] <- 1

hyp2a <- lm(CuresMisinfo ~ treatment * religiosity, data=hyp2)

hyp2b <- lm(ConspiracyMisinfo ~ treatment  * religiosity, data=hyp2)

stargazer(hyp2b, hyp2a, star.cutoffs =  c(0.05, .001, 0.001), type="text")

##### table G2 #####
hyp3 <- subset(covid, covid$B5_Q1==2 | covid$B5_Q1==5)
hyp3$treatment <- 0
hyp3$treatment[hyp3$B5_Q1==2] <- 1

hyp3a <- lm(CuresMisinfo ~ treatment  * religiosity, data=hyp3)

hyp3b <- lm(ConspiracyMisinfo ~ treatment  * religiosity, data=hyp3)

stargazer(hyp3b, hyp3a, star.cutoffs =  c(0.05, .001, 0.001), type="text")

##### table H1 and H2 #####

# table H1
hyp1b <- lm(AllMisinfo ~ religiosity, data=covid)
stargazer(hyp1b, star.cutoffs =  c(0.05, .001, 0.001), type="text")

# table H2
hyp1b_covars <- lm(AllMisinfo ~ religiosity + BJP +
                     B1_Q2 + gender + B1_Q8 + education + uppercaste +
                     science_knowledge + B2_Q7 + B2_Q8, data=covid)

stargazer(hyp1b_covars, star.cutoffs =  c(0.05, .001, 0.001), type="text")


##### table I1 #####

# conspiracies: pooled dv minus the corrected story (sneezing)
covid$B5_Q1 <- factor(covid$B5_Q1, levels = c("5", "1", "2", "3", "4"))

covid$ConspiracyRobust <- covid$jamaat + covid$bio + covid$foreign +
  covid$lab + covid$fiveg
robust_cons <- lm(ConspiracyRobust ~ factor(B5_Q1), data=covid)

# cures: pooled dv minus the corrected story (ayurveda)
covid$CuresRobust <- covid$kalonji + covid$genetic + covid$breath +
  covid$antib + covid$bleach
robust_cures <- lm(CuresRobust ~ factor(B5_Q1), data=covid)

stargazer(robust_cons, robust_cures, star.cutoffs = c(0.05, .001, 0.001), type="text")


##### table I2 #####

# conspiracies: pooled dv minus the corrected story (sneezing)
covid$B5_Q1 <- factor(covid$B5_Q1, levels = c("4", "1", "2", "3", "5"))

covid$ConspiracyRobust2 <- covid$jamaat + covid$bio + covid$foreign +
  covid$lab + covid$fiveg
robust_cons2 <- lm(ConspiracyRobust2 ~ factor(B5_Q1), data=covid)

# cures: pooled dv minus the corrected story (ayurveda)
covid$CuresRobust2 <- covid$kalonji + covid$genetic + covid$breath +
  covid$antib + covid$bleach
robust_cures2 <- lm(CuresRobust2 ~ factor(B5_Q1), data=covid)

stargazer(robust_cons2, robust_cures2, star.cutoffs = c(0.05, .001, 0.001), type="text")


##### table J1 #####
pol_all1 <- lm(AllMisinfo ~ B3_Q12, data=covid)
stargazer(pol_all1, star.cutoffs =  c(0.05, .001, 0.001), type="text")


##### figure K2 #####

# astrology
covid$relig5 <- as.factor(covid$relig5)
levels(covid$relig5) <- c("Strongly Disagree", "Disagree", "Agree", "Strongly Agree")
astrology <- lm(AllMisinfo ~ relig5, data=covid)
summary(astrology)

#vegetarian food
covid$relig6 <- as.factor(covid$relig6)
levels(covid$relig6) <- c("Strongly Disagree", "Disagree", "Agree", "Strongly Agree")
vegfood <- lm(AllMisinfo ~ relig6, data=covid)
summary(vegfood)

#puja god blesses me
covid$relig4 <- as.factor(covid$relig4)
levels(covid$relig4) <- c("Strongly Disagree", "Disagree", "Agree", "Strongly Agree")
puja <- lm(AllMisinfo ~ relig4, data=covid)
summary(puja)

#fasting
covid$relig3 <- as.factor(covid$relig3)
levels(covid$relig3) <- c("Strongly Disagree", "Disagree", "Agree", "Strongly Agree")
fasting <- lm(AllMisinfo ~ relig3, data=covid)
summary(fasting)

#important to teach children
covid$relig2 <- as.factor(covid$relig2)
levels(covid$relig2) <- c("Strongly Disagree", "Disagree", "Agree", "Strongly Agree")
teach <- lm(AllMisinfo ~ relig2, data=covid)
summary(teach)

#religion helps me cope
covid$relig1 <- as.factor(covid$relig1)
levels(covid$relig1) <- c("Strongly Disagree", "Disagree", "Agree", "Strongly Agree")
cope <- lm(AllMisinfo ~ relig1, data=covid)
summary(cope)

#athiest can be moral (was reverse coded)
covid$B3_Q8 <- as.factor(covid$B3_Q8)
levels(covid$B3_Q8) <- c("Strongly Agree", "Agree", "Disagree", "Strongly Disagree")
athiest <- lm(AllMisinfo ~ B3_Q8, data=covid)
summary(athiest)

#marry someone not hindu 
covid$B3_Q1 <- as.factor(covid$B3_Q1)
levels(covid$B3_Q1) <- c("Strongly Agree", "Agree", "Disagree", "Strongly Disagree")
table(covid$B3_Q1)
marry <- lm(AllMisinfo ~ B3_Q1, data=covid)
summary(marry)

# run below line if R tells you margins are off

#par(mar=c(2.5,2.5,2.5,2.5)) 

par(mfrow=c(4,2))

cplot(astrology, "relig5", what = "prediction", 
      main = "Consult astrologer before wedding", 
      xlab="", 
      ylab="# of stories accurately identified",
      ylim=c(4.5,8))

cplot(vegfood, "relig6", what = "prediction", 
      main = "Eat vegetarian food only", 
      xlab="", 
      ylab="# of stories accurately identified",
      ylim=c(4.5,8))

cplot(puja, "relig4", what = "prediction", 
      main = "God blesses me when I pray", 
      xlab="", 
      ylab="# of stories accurately identified",
      ylim=c(4.5,8))

cplot(fasting, "relig3", what = "prediction", 
      main = "Fasting is important", 
      xlab="", 
      ylab="# of stories accurately identified",
      ylim=c(4.5,8))


cplot(teach, "relig2", what = "prediction", 
      main = "Imp. to teach children religion", 
      xlab="", 
      ylab="# of stories accurately identified",
      ylim=c(4.5,8))

cplot(cope, "relig1", what = "prediction", 
      main = "Religion helps me cope", 
      xlab="", 
      ylab="# of stories accurately identified",
      ylim=c(4.5,8))

cplot(athiest, "B3_Q8", what = "prediction", 
      main = "Athiests can be moral", 
      xlab="", 
      ylab="# of stories accurately identified",
      ylim=c(4.5,8))

cplot(marry, "B3_Q1", what = "prediction", 
      main = "Marry someone not Hindu", 
      xlab="", 
      ylab="# of stories accurately identified",
      ylim=c(4.5,8))


##### table L1 and L2 #####
covid$B5_Q1 <- factor(covid$B5_Q1, levels = c("5", "1", "2", "3", "4"))
levels(covid$B5_Q1)

# conspiracy stories

cons1 <- lm(jamaat ~ factor(B5_Q1), data=covid)
summary(cons1)
cons2 <- lm(sneezing ~ factor(B5_Q1), data=covid)
summary(cons2)
cons3 <- lm(bio ~ factor(B5_Q1), data=covid)
summary(cons3)
cons4 <- lm(foreign ~ factor(B5_Q1), data=covid)
summary(cons4)
cons5 <- lm(lab ~ factor(B5_Q1), data=covid)
summary(cons5)
cons6 <- lm(fiveg ~ factor(B5_Q1), data=covid)
summary(cons6)

stargazer(cons1, cons2, cons3, cons4, cons5, cons6, 
          star.cutoffs =  c(0.05, .001, 0.001), type="text")

# medical stories

cure1 <- lm(homeo ~ factor(B5_Q1), data=covid)
summary(cure1)
cure2 <- lm(kalonji ~ factor(B5_Q1), data=covid)
summary(cure2)
cure3 <- lm(genetic ~ factor(B5_Q1), data=covid)
summary(cure3)
cure4 <- lm(breath ~ factor(B5_Q1), data=covid)
summary(cure4)
cure5 <- lm(antib ~ factor(B5_Q1), data=covid)
summary(cure5)
cure6 <- lm(bleach ~ factor(B5_Q1), data=covid)
summary(cure6)

stargazer(cure1, cure2, cure3, cure4, cure5, cure6, 
          star.cutoffs =  c(0.05, .001, 0.001), type="text")

##### table M1 #####

covid$B5_Q1 <- factor(covid$B5_Q1, levels = c("1", "2", "3", "4", "5"))

main_cures <- lm(CuresMisinfo ~ factor(B5_Q1), data=covid)

main_cons <- lm(ConspiracyMisinfo ~ factor(B5_Q1), data=covid)

stargazer(main_cons, main_cures, star.cutoffs =  c(0.05, .001, 0.001), type="text")

##### table N1 #####

covid$religion_themed <- covid$jamaat + covid$sneezing

covid$not_religion_themed <- covid$bio + covid$foreign + covid$lab + covid$fiveg

covid$B5_Q1 <- factor(covid$B5_Q1, levels = c("5", "1", "2", "3", "4"))

newreg1 <- lm(religion_themed ~ factor(B5_Q1), data=covid)

newreg2 <- lm(not_religion_themed ~ factor(B5_Q1), data=covid)

stargazer(newreg1, newreg2, star.cutoffs =  c(0.05, .001, 0.001), type="text")

##### table O1 ##### 
hyp1a <- lm(AllMisinfo ~ BJP,  data=covid)
summary(hyp1a)
stargazer(hyp1a, star.cutoffs =  c(0.05, .001, 0.001), type="text")

##### table O2 #####

table(covid$B5_Q1)

hyp4 <- subset(covid, covid$B5_Q1==3 | covid$B5_Q1==5)
table(hyp4$B5_Q1)
hyp4$treatment <- 0
hyp4$treatment[hyp4$B5_Q1==3] <- 1
table(hyp4$treatment)

hyp4a <- lm(CuresMisinfo ~ treatment * BJP, data=hyp4)
summary(hyp4a)

hyp4b <- lm(ConspiracyMisinfo ~ treatment *  BJP, data=hyp4)
summary(hyp4b)

stargazer(hyp4b, hyp4a, star.cutoffs =  c(0.05, .001, 0.001), type="text")



